Real time radiation treatment planning system

ABSTRACT

The invention relates to a real time radiation treatment planning system for use in effecting radiation therapy of a pre-selected anatomical portion of an animal body using hollow needles. According to embodiments of the invention, the system may include a processing means processing means-configured to perform a three-dimensional imaging algorithm and a three-dimensional image segmentation algorithm, with respect to one or more specific organs within the pre-selected anatomical portion and with respect to the needles, for converting the image data obtained with an imaging means into a three-dimensional image of the anatomical portion, using at least one single or multi-objective anatomy-based genetic optimization algorithm. For pre-planning or virtual simulation purposes, the processing means is arranged to determine in real time the optimal number and position of at least one of the needles, positions of energy emitting sources within the needles, and the dwell times of the energy emitting sources at the positions. For post-planning purposes, the processing means is arranged to determine, based on three-dimensional image information, in real time the real needle positions and the dwell times of the energy emitting sources for the positions.

This application is a continuation application of U.S. application Ser.No. 10/448,311, entitled REAL TIME RADIATION TREATMENT PLANNING SYSTEM,filed May 30, 2003, of Johann KINDLEIN et al., currently pending, whichclaims the right of priority under 35 U.S.C. §119 to EuropeanApplication No. 02077436.0 filed on Jun. 17, 2002, the entiredisclosures of which are fully incorporated herein by reference.

BACKGROUND

The invention relates to a real time radiation treatment planning systemfor use in effecting radiation therapy of a pre-selected anatomicalportion of an animal body, comprising:

A a stepper for automatically positioning imaging means for generatingimage data corresponding to the anatomical portion;B means for inserting under the guidance of a template at least onehollow needle at a position into said anatomical portion;C radiation delivery means for defining a plurality of positions havinga spatial relationship within a volume of said anatomical portion andfor inserting at least one energy emitting source through said at leastone hollow needle at said plurality of positions into said anatomicalportion;D processing means for generating a radiation treatment plan foreffecting said radiation therapy, said treatment plan includinginformation concerning:

the number, position, direction and estimation of the best way ofplacement of one or more of said hollow needles within the anatomicalportion and volume of said anatomical portion to be treated;

the amount of radiation dose to be emitted.

The last decade has seen major changes in the way radiation treatmentsare delivered. The century-old objective of radiation therapy, i.e. todeliver a curative dose to the target, e.g. a tumour, while preservingnormal tissues of the animal body can now be aimed at with a high degreeof sophistication. However, despite of major improvements achieved withthree-dimensional imaging techniques, that allow the anatomy to beproperly defined, brachytherapy treatments have not yet fully benefitedfrom these important new pieces of information.

For brachytherapy using high dose rate (HDR) energy emitting sources,catheters or hollow needles are placed in a target volume within ananimal body and it is assumed that if the dose distribution covers thecatheters, it should also cover the anatomy. Imaging is commonly used toset the treatment margins, but optimized dose distributions are based onconsiderations, such as the catheter positions and desired dose andlimited to a few defined points. This necessarily results in anapproximation of the shape of the anatomical portion to be treated.

For the case of treatments of the prostate, volume optimization resultsin a dose distribution that is essentially cylindrically shaped. With acylindrically shaped approximation of the prostate it is possible toassure the complete coverage of the prostate volume with the radiationemitted by the source or sources. Only a conformal dose distributiondelivered to the anatomical portion with an adequate margin around theprostate will encompass all affected, cancerous tissue.

The methods described in the prior art (e.g. Etienne Lessard, Med. Phys.28. (5), May 2001) are using the concept of inverse planning to obtainan anatomy-based optimization of the dose distribution. Without anymanual modification to deliver conformal HDR prostate treatment andknowing the exact location of the applicators (catheters/hollowneedles), due to modern imaging techniques, it is easy to determine thepossible stopping position of the radioactive source within a catheteror hollow needle present in the animal body. The possible sourcepositions are considered given. The system has to determine based on aHDR inverse planning dose optimization governed entirely from anatomyand clinical criteria to decide the best dwell time distribution.

In U.S. Pat. No. 5,391,139 in the name of G. K. Edmundson a real timeradiation treatment planning system according to the preamble isdisclosed. With this system image data of the anatomical portion, e.g.the prostate is obtained for planning purposes and the medical personnelchooses an arbitrary number of needle locations using predeterminedplacement rules, which have been empirically determined from experience.The planning system develops a treatment plan based on these arbitraryneedle positions after which the medical personnel has to examine theplanning results and decide whether these results are suitable for theperforming the actual radiation treatment. In case the medical personnelfinds the planning results unsatisfactorily the virtual needle positionshave to be altered and using the repositioned needles a new treatmentplan is generated. This trial-and-error approach is repeated until atreatment plan is developed that satisfies the actual intended radiationtreatment.

Subsequently the catheters or needles are inserted via a template intothe animal body according to the generated treatment plan.

Conventional dose optimization algorithms are single objective, i.e.they provide a single solution. This solution is found by atrial-and-error search method as in Edmundson's U.S. Pat. No. 5,391,139,by modifying importance factors of a weighted sum of objectives, e.g. byrepositioning the virtual needles or by changing the radiation dose tobe delivered. This problem has been addressed currently and some methodshave been proposed to find an optimal set of importance factors.

Conventional optimization methods combine the target objectives and theobjectives for the surrounding healthy tissue and of critical structuresinto a single weighted objective function. The weight or importancefactor for each objective must be supplied. The obtained solutiondepends on the value of importance factors used. One goal of a treatmentplanning system is the ability to assist the clinician in obtaining goodplans on the fly. Also it should provide all the information of thepossibilities given the objectives of the treatment. In order to explorethe feasible region of the solution space with respect to eachobjective, different values for the importance factors in the aggregateobjective function must be given.

Furthermore, the appropriate values of these importance factors differfrom clinical case to clinical case. This implies that for any newclinical case a lot of effort is necessary for their determination.

While current optimization methods are single weighted objective methodsthe dose optimization problem is a true multi-objective problem andtherefore multi-objective optimization methods should be used.

The gradient-based algorithm due to its efficiency allows theconstruction of the so-called Pareto or trade-off surface which containsall the information of the competition between the objectives which isnecessary for the planner to select the solution which best fulfills hisrequirements.

One problem of this algorithm is that the weighted sum as used in allconventional dose optimization algorithms cannot provide solutions inpossible non-convex parts of the Pareto tradeoff surface, because aconvex weighted sum of objectives converges only to the convex parts ofthe Pareto front. Another major limitation of the algorithm is itsrestriction to convex objective functions for which gradients can becalculated. In this case according to the Kuhn-Tucker theorems a globaloptimum can be obtained and the entire Pareto front is accessible fromthe weighted sum.

When searching for an optimal set of importance factors dividing eachimportance factors in n points, then the number of combinations for kobjectives is approximately proportional to n^(k-1) and the shape of theentire trade-off surface require a very large computational time. Mostrealistic problems require the simultaneous optimization of manyobjectives. It is unlikely that all objectives are optimal for a singleset of parameters. If this is so, then there exist many, in principleinfinite solutions.

A multi-objective algorithm does not provide a single solution, but arepresentative set of all possible solutions. Out of theserepresentative solutions a single final solution has to be selected. Itis a complex problem to automatically select such a solution and suchmethods have been proposed but then a planner would not know whatalternatives solutions could instead be selected. In problems wheredifferent sets of objectives have to be compared this information isvaluable, since it shows the possibilities a planner has for each suchset.

A time analysis of the optimization with available commercial systemsbased on e.g. 35 clinical cases shows that even if a single optimizationrun requires only a few seconds the actual optimization requires 5.7±4.8minutes. The evaluation of the results requires additional 5.8±2.5minutes. This shows that the result of a single optimization run is notalways satisfactorily and most of the time is spent in a manualtrial-and-error optimization.

The invention aims to obviate the above described problems and proposesa new real time radiation treatment planning system according to theabove preamble, where the possible positions of the energy emittingsources are not considered given and the location of the needles are notpredetermined based on rules, which have been empirically determined.

It is also an object of the present invention to describe a new realtime radiation treatment planning system and method, which will allow asignificant speed-up of single and multi-objective anatomy based doseoptimization and inverse planning procedures for HDR interstitialbrachytherapy.

More in particular the invention aims to generate in real time atreatment plan, which will be presented to the medical personnelinstantly and can also immediately used as the radiation treatment.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 shows a schematic of an exemplary device for implanting an energyemitting source;

FIG. 2 shows an exemplary template for use in a real time radiationtreatment planning system;

FIG. 3 shows an exemplary catheter or hollow needle;

FIG. 4 shows exemplary generated sampling points distributed on thecontours of the anatomical portion to be treated;

FIG. 5 shows exemplary generated sampling points distributed on thetriangulated surface of the anatomical portion to be treated;

FIG. 6 shows an exemplary dialog for displaying the template andcatheter characteristics;

FIG. 7 shows an exemplary “auto-activation” dialog;

FIG. 8 shows an exemplary “source parameter” dialog;

FIG. 9 shows an exemplary “inverse planning” dialog;

FIG. 10 shows an exemplary “template view and loading” dialog;

FIG. 11 shows three exemplary “loading” buttons, consistent with theembodiment of FIG. 10;

FIG. 12 shows an exemplary “geometry and sampling” dialog;

FIG. 13 shows an exemplary “optimization” dialog;

FIG. 14 shows an exemplary “post-plan dose optimization” dialog;

FIG. 15 shows an exemplary “analysis of solutions” dialog;

FIG. 16 shows a first exemplary view of a “show results of allsolutions” dialog;

FIG. 17 shows a second exemplary of the “show results of all solutions”dialog;

FIG. 18 shows a third exemplary view of the “show results of allsolutions” dialog with values in DVH column sorted in descending order;

FIG. 19 shows exemplary histograms of radiation dose coverage of theanatomical portion (PTV);

FIG. 20 shows an exemplary “analysis of solutions” dialog;

FIG. 21 shows an exemplary “dose distribution” dialog;

FIG. 22 shows another exemplary view of the “auto-activation” dialog;and

FIG. 23 shows an exemplary “source parameter” dialog.

DETAILED DESCRIPTION

According to the invention said processing means are provided with athree-dimensional imaging algorithm and a three-dimensional imagesegmentation algorithm for at least the specific organs within saidanatomical portion and the needles for converting the image dataobtained with said imaging means into a three-dimensional image of theanatomical portion, whereby by using at least one single ormulti-objective anatomy based genetic optimization algorithm forpre-planning or virtual simulation purposes said processing means arearranged to determine in real time the optimal number and position of atleast one of said hollow needles, the position of said energy emittingsource within each hollow needle as well as the dwell times of saidenergy emitting source at each position; whereas for post planningpurposes said processing means are arranged to determine based onthree-dimensional image information in real time the real needlepositions and the dwell times of said energy emitting source for eachposition.

The term needles in this application also covers e.g. catheters, guidetubes or other elements to be implanted in an animal body forpositioning an energy emitting source inside that animal body. Thethree-dimensional image segmentation algorithm used in the treatmentplanning system according to the invention may use all or part of theseelements (hollow needles, catheters, and guide tubes, e.g.).

The use of a three-dimensional imaging algorithm and a three-dimensionalimage segmentation algorithm allows a real time and fast generation of atreatment plan without the necessity of determining certain objectivesas a starting point for the treatment planning step, such as thepositioning of one or more needle inside the anatomical portion.

In fact with the treatment planning system according to the inventionfor generating a treatment plan said anatomy based genetic optimizationalgorithm uses specific animal body related data and/or system relateddata and/or radiation related data, wherein said animal body relateddata are data concerning the shape and localization of said anatomicalportion and/or the shape and localization of specific organs near orwithin said anatomical portion.

Said system related data may be data concerning the template and itsposition in relation to said anatomical portion of said animal bodyand/or the dimensions of the needles used and/or the minimumdisplacement distance of said energy emitting source through saidradiation delivery means, whereas said radiation related data are dataconcerning the prescribed radiation dose to said anatomical portion ofsaid animal body, the maximum radiation exposure dose to said specificorgans near or within said anatomical portion.

These specific data are used as boundary conditions and entered into theanatomy based genetic optimization algorithm by the medical personnel ordetermined/established by said algorithm from e.g. said image dataobtained with the imaging means.

Dose optimization has to consider many objectives in conflict, such asthe coverage of the pre-selected anatomical portion or planning targetvolume (PTV) to be treated with a specified dose and the dose protectionof the surrounding tissue and specific delicate organs (OAR=organs atrisk), such as the bladder and urethra, when treating prostate cancer.The objectives are combined into a single objective function ƒ_(Tot)formed by a weighted sum of the individual objective functions. Theoptimal value f*i for the i^(th) objective found by a optimizationalgorithm depends on the weights (importance factors) used and may notbe the best possible result as the mapping from importance to objectivespace is complex, especially for three and more objectives. In caseswhere the solution is not satisfactory the treatment planner is requiredto repeat the optimization with a different set of importance factors.One method is to increase the importance factors of the objectives forwhich the solution does not provide a satisfactory result.

Practical only a very small number of combinations can be tested andwith this approach the treatment planner cannot gain all the informationabout the range of possible values and the degree of competition, whichare required to select the “best solution”. In order to get the “best”possible result avoiding trial-and-error methods the invention proposesand uses a gradient based multi-objective optimization algorithm.

The importance of this has been recognized by the U.S. Pat. No.6,327,490, where a method is presented to save and compare differenttreatment plans. What is not recognized in this patent is that not onlya modification of the energy emitting source positions is sometimesnecessary, but also a modification of the importance factors/boundaryconditions, which determine the quality of the solution. In clinicalpractice the reality is often that the optimization required ismulti-objective rather than single-objective and this increases evenfurther the time required for the calculation process. Consequentlythere is a real need in anatomy based optimization techniques andinverse planning procedures for very fast dose calculation methods.

The dose d_(i)(X) at the i^(th) sampling point is calculated by:

${d_{i}(x)} = {\sum\limits_{j = 1}^{N_{d}}{x_{j}^{2}{\overset{\sim}{d}}_{{ij}\;}}}$

where N_(d) is the number of sources, x_(j) ² is the dwell time of thej^(th) source dwell position and the kernel value for the i^(th) dosecalculation point and j^(th) source dwell position. Dose calculationlook-up tables (LUT) of {tilde over (d)}_(ij) are calculated and storedin a preprocessing step. The calculation of the dose for N_(s) samplingpoints requires N_(s)N_(d) multiplications and N_(s)(N_(d)−1) additions.

A Fast Fourier Transform (FFT) based convolution technique for thecalculation of a three-dimensional dose distribution has been proposedby Boyer and Mok in 1986, which tries to reduce the calculation time.Employing FFT based convolution methods the time needed for thecalculation of a dose distribution with this method is independent ofthe number of sources. While unable to cope with angular asymmetrickernels, an advantage of the FFT based method is that the computationaltime is practically independent of the form of the dosimetric kernelused. Our analysis has shown that this method is comparable withconventional methods only if the number of sources is much larger than300. The avoidance of wrap around effects can only be avoided with zeropadding in each dimension which increases the transform size toN=8.N_(s) for N_(s) sampling points and requires 8.N_(s).ln(8.N_(s))operations.

The aim of HDR brachytherapy dose optimization according to theinvention is to cover the anatomical portion (PTV) to be treated with atleast some dose value, and to protect the specific delicate organs (OAR)and the surrounding normal tissue with dose values above some specificlevel. For variance based objective functions dose values above acritical dose value are penalized quadratic. The objectives are suchthat the isodose of the optimal dose distribution of the prescriptiondose coincides with the surface of the anatomical portion. With thisapproach, the use of an additional objective for the surrounding normaltissue is not necessary. For this the dose variance f_(s) of thesampling points (dose points) as uniformly distributed on the surface ofthe anatomical portion should be as small as possible. The avoidance ofexcessive high dose values inside the anatomical portion, e.g. theprostate is controlled by the dose distribution variance f_(v) insidethe anatomical portion.

Normalized variances are used:

${f_{s} = {\frac{1}{N_{s}}{\sum\limits_{i = 1}^{N_{s}}\frac{\left( {d_{i\;}^{s} - m_{s}} \right)^{2}}{m_{s}^{2}}}}},{f_{v} = {\frac{1}{N_{v}}{\sum\limits_{i = 1}^{N_{v\;}}\frac{\left( {d_{i}^{v} - m_{v}} \right)^{2}}{m_{v}^{2}}}}}$

where m_(s) and m_(v) are the average dose values on the surface of theanatomical portion and within the anatomical portion respectively, andN_(s), N_(v) the corresponding numbers of sampling points. The objectivespace of (f_(s), f_(v)) is convex and gradient-based algorithms convergeto the global Pareto front. If specific, delicate organs are to beconsidered then an additional objective is included for each specificorgan (OAR):

${f_{OAR} = {\frac{1}{N_{OAR}}{\sum\limits_{i = 1}^{N_{OAR}}\frac{{\Theta \left( {d_{i}^{OAR} - {D_{c}^{OAR}m_{s}}} \right)}\left( {d_{i\;}^{OAR} - {D_{c}^{OAR}m_{s}}} \right)^{2}}{\left( {D_{c}^{OAR}m_{s}} \right)^{2}}}}},{{\Theta (x)} = \left\{ \begin{matrix}1 & {x \geq 0} \\0 & {x < 0}\end{matrix} \right.}$

where N_(OAR) is the number of sampling points in the specific organ andD_(c) ^(OAR) is the corresponding critical dose as a fraction of theprescription dose or reference dose, which dose equals in this model theaverage dose on the surface of the anatomical portion. The objectivefunctions for the specific delicate organs are of the same form as forthe anatomical organs, but involve the dose variances versus thecritical dose values, which are specific only to those particularspecific organs. The derivatives are:

$\frac{\partial f_{s\;}}{\partial x_{k}} = {\frac{4x_{k}}{N_{s}m_{s}^{3}}{\sum\limits_{i = 1}^{N_{s}}{d_{i}^{s}\left( {{m_{s}{\overset{\sim}{d}}_{ik}^{s}} - {d_{i}^{s}{\overset{\sim}{m}}_{k}^{s}}} \right)}}}$$\frac{\partial f_{v}}{\partial x_{k}} = {\frac{4x_{k}}{N_{v}m_{v\;}^{3}}{\sum\limits_{i = 1}^{N_{v\;}}{d_{i}^{v}\left( {{m_{v}{\overset{\sim}{d}}_{ik}^{v}} - {d_{i}^{v}{\overset{\sim}{m}}_{k}^{v}}} \right)}}}$$\frac{\partial f_{OAR}}{\partial x_{k}} = {\frac{4x_{k}}{{N_{OAR}\left( D_{C}^{OAR} \right)}^{2}m_{s}^{3}}{\sum\limits_{i = 1}^{N_{OAR}}{{\Theta \left( {d_{i}^{OAR} - {D_{c}^{OAR}m_{s}}} \right)}\left( {d_{i}^{OAR} - {D_{c}^{OAR}m_{s}}} \right)\left( {{m_{s}{\overset{\sim}{d}}_{ik}^{OAR}} - {d_{i}^{OAR}{\overset{\sim}{m}}_{k}^{s}}} \right)}}}$

Where the following relations are used:

${d_{i}^{s} = {\sum\limits_{i = 1}^{N_{d}}{x_{i}^{2}{\overset{\sim}{d}}_{il}^{s}}}},{m_{s} = {\frac{1}{N_{s}}{\sum\limits_{i = 1}^{N_{s}}d_{l}^{s}}}},{{\overset{\sim}{m}}_{k}^{s} = {\frac{1}{N_{s}}{\sum\limits_{i = 1}^{N_{s}}{\overset{\sim}{d}}_{lk}^{s}}}},{k = 1},{\ldots \mspace{14mu} N_{d}}$${d_{i}^{v} = {\sum\limits_{i = 1}^{N_{d}}{x_{l}^{2}{\overset{\sim}{d}}_{il}^{v}}}},{m_{v} = {\frac{1}{N_{v}}{\sum\limits_{i = 1}^{N_{v}}d_{l}^{v}}}},{{\overset{\sim}{m}}_{k}^{v} = {\frac{1}{N_{v}}{\sum\limits_{i = 1}^{N_{v}}{\overset{\sim}{d}}_{lk}^{v}}}},{k = 1},{\ldots \mspace{14mu} N_{d}}$$d_{i}^{OAR} = {\sum\limits_{i = 1}^{N_{d}}{x_{l}^{2}{\overset{\sim}{d}}_{il}^{OAR}}}$

where d_(i) ^(s), d_(i) ^(v) and d_(i) ^(OAR) is the dose rate at thei^(th) sampling point on the surface of the anatomical portion, withinthe anatomical portion and within a specific organ respectively. {tildeover (d)}_(i) _(il) ^(S), {tilde over (d)}_(i) _(il) ^(V), {tilde over(d)}_(i) _(il) ^(OAR) is the dose kernel for the i^(th) sampling pointand the l^(th) source dwell position for the sampling points on thesurface of the anatomical portion, within the anatomical portion and inthe specific organ respectively. N_(d) is the number of source dwellpositions.

For a conventional method to calculate objective values and derivativesthe lookup table {tilde over (d)}_(i) _(ia) ^(S) as a N_(s)×N_(d) matrix{tilde over (K)}^(S) are considered:

$K^{s} = \begin{pmatrix}{\overset{\sim}{d}}_{il}^{s} & {\overset{\sim}{d}}_{12}^{s} & \ldots & {\overset{\sim}{d}}_{1{Nd}}^{s} \\{\overset{\sim}{d}}_{21}^{s} & {\overset{\sim}{d}}_{21}^{s} & \ldots & {\overset{\sim}{d}}_{2{Nd}}^{s} \\\vdots & \vdots & \ddots & \vdots \\{\overset{\sim}{d}}_{Nsl}^{s} & {\overset{\sim}{d}}_{{Ns}\; 2}^{s} & \ldots & {\overset{\sim}{d}}_{NsNd}^{s}\end{pmatrix}$

If d^(s) ^(T) =(d₁ ^(s), d₂ ^(s), . . . , d_(Ns) ^(s)) is the vector ofdose values d_(i) ^(s) and t^(−T)=(t₁, t₂, . . . , t_(Nd))=(x₁ ², x₂ ²,. . . , x_(Nd) ²) the vector of the dwell times, then d^(s)={tilde over(K)}^(s) t.

The conventional approach is to calculate the dose values and then theobjective values and their derivatives. The number of operations tocalculate d ^(s) requires N_(d).N_(S) multiplications and(N_(d)−1).N_(s) additions. The storage required for {tilde over (K)}^(S)is N_(d).N_(s) floating points. The storage of this pre-computed matrixis desired because of the significant gain in the optimization speed.The calculation of the objective functions and the derivatives for N_(s)sampling points and N_(d) source dwell positions requires therefore anorder of N_(d).N_(s) operations.

The New Method of the Calculation of the Objectives and DerivativesUsing Dose Kernel Look-Up Tables

The objective function ƒ_(s) and its N_(d) derivatives can be writtenalso as:

${f_{s} = {\left\lbrack {\frac{1}{N_{s}m_{s}^{2}}{\sum\limits_{i = l}^{N_{s}}\left( d_{1}^{s} \right)^{2}}} \right\rbrack - 1}},{\frac{\partial f_{s}}{\partial x_{k}} = {{{\frac{4x_{k}}{N_{s}m_{s}^{2}}{\sum\limits_{i = 1}^{N_{s}}{d_{i}^{s}{\overset{\sim}{d}}_{ik}^{s}}}} + {\frac{4x_{k}{\overset{\sim}{m}}_{k}^{s}}{N_{s}m_{s}^{3}}{\sum\limits_{i = 1}^{N_{s}}{\left( d_{i}^{s} \right)^{2}k}}}} = 1}},2,\ldots \mspace{14mu},N_{d}$

It is possible to reduce the number of operations of the objective valueand its derivatives to approximately O(N_(d).N_(d)) operations. With{right arrow over (d)}^(S) ^(T) ={right arrow over (t)}^(T) {tilde over(K)}^(S) ^(T) we have

${\sum\limits_{i = 1}^{NS}\left( d_{i}^{S} \right)^{2}} = {{{\overset{\rightarrow}{d}}^{S^{T}} \cdot {\overset{\rightarrow}{d}}^{S}} = {{{\overset{\rightarrow}{t}}^{T}{\overset{\sim}{K}}^{S^{T}}{\overset{\sim}{K}}^{S}\overset{\rightarrow}{t}} = {{\overset{\rightarrow}{t}}^{T}D^{s}{\overset{\rightarrow}{t}.}}}}$

D^(s) is a symmetric (D_(αβ) ^(S)=D_(βα) ^(S)) N_(d)×N_(d) matrix. Wecall {tilde over (K)}^(S) the first order and D^(s) the second orderdose kernel matrix.

The terms m_(s), m_(s) ² and

$\sum\limits_{i = 1}^{N_{S}}{d_{i}^{S}{\overset{\sim}{d}}_{i\; \beta}^{S}}$

can be calculated from {tilde over (K)}^(S) using the relations:

${m_{S} = {{\overset{\rightarrow}{s}}^{S^{T}}{\overset{\sim}{K}}^{S}\overset{\rightarrow}{t}}},{m_{S}^{2} = {{{\overset{\rightarrow}{t}}^{T}{\overset{\sim}{K}}^{S^{T}}{\overset{\rightarrow}{s}}^{S}{\overset{\rightarrow}{s}}^{S^{T}}{\overset{\sim}{K}}^{S}\overset{\rightarrow}{t}} = {{\overset{\rightarrow}{t}}^{T}{\overset{\overset{\rightarrow}{\sim}}{m}}^{S}{\overset{\overset{\rightarrow}{\sim}}{m}}^{S^{T}}\overset{\rightarrow}{t}}}},{{\sum\limits_{i = 1}^{N_{S}}{d_{i}^{S}{\overset{\sim}{d}}_{i\; \beta}^{S}}} = \left( {{\overset{\rightarrow}{s}}^{S^{T}}{\overset{\sim}{K}}^{S}} \right)_{\beta}}$

where

${\overset{\rightarrow}{s}}^{S^{T}} = {\frac{1}{N_{S}}\left( {1,1,1,\ldots \mspace{14mu},1} \right)}$

is a N_(s) dimensional vector and {tilde over ({right arrow over(m)}^(S) ^(T) ={right arrow over (s)}^(S) ^(T) {tilde over (K)}^(S).

From the matrix representation the following terms can be writtenanalytically as:

$\begin{matrix}{{g_{\alpha} \equiv {\sum\limits_{i = 1}^{N_{S}}{d_{i}^{S}{\overset{\sim}{d}}_{i\; \alpha}^{S}}}} = {\sum\limits_{\beta = 1}^{N_{d}}{x_{\beta}^{2}D_{\alpha\beta}^{S}}}} & (1) \\{{\sum\limits_{i = 1}^{N_{S}}\left( d_{i}^{S} \right)^{2}} = {{\sum\limits_{\alpha = 1}^{N_{d}}{\sum\limits_{\beta = 1}^{N_{d}}{x_{\alpha}^{2}x_{\beta}^{2}D_{\alpha \; \beta}^{S}}}} = {{\sum\limits_{\alpha = 1}^{N_{d}}{x_{\alpha}^{2}{\sum\limits_{\beta = 1}^{N_{d}}{x_{\beta}^{2}D_{\alpha \; \beta}^{S}}}}} = {\sum\limits_{\alpha = 1}^{N_{d}}{x_{\alpha}^{2}g_{\alpha \;}}}}}} & (2) \\{m_{S} = {\sum\limits_{\alpha = 1}^{N_{d}}{x_{\alpha}^{2}{\overset{\sim}{m}}_{\alpha}^{S}}}} & (3)\end{matrix}$

where

$D_{\alpha \; \beta}^{S} = {\sum\limits_{i = 1}^{N_{S}}{{\overset{\sim}{d}}_{i\; \alpha}^{S}{{\overset{\sim}{d}}_{i\; \beta}^{S}.}}}$

For the objective f_(v) the corresponding matrices {tilde over (K)}^(V)and {tilde over (D)}^(v) are required. D^(v) and D^(v) can be calculatedonce in a preprocessing step which requires only 1-2 seconds.

The equations (1)-(3) show us that it is possible to increase the numberof sampling points in the anatomical portion without increasing theoptimization time as the right side of each equation does not directlydepend on the number of sampling points. This means that we can increaseaccuracy without increasing the computation costs.

With the new approach the number of operations is independent on thenumber of sampling points. It is not necessary to store the matrix{tilde over (K)}^(S) and {tilde over (K)}^(V). Only the matrices D^(S),D^(V) and the N_(d) dimensional vectors {tilde over (K)}^(sT) s ^(s) and{tilde over (K)}^(vT) s ^(v) are required for the dose optimization. Thematrices {tilde over (K)}^(S) and {tilde over (K)}^(V) requireN_(s).N_(d) and N_(v).N_(d) numbers to be stored, whereas both symmetricmatrices D^(S) and D^(V) require N_(d) (N_(d)+1) numbers to be stored,i.e. the storage is independent on the number of sampling points.

It is of importance to mention that this novel algorithm can be appliedin general for objectives (and their derivatives) of the following typescommonly used in brachytherapy (HDR and seeds) and external beamradiotherapy:

$\begin{matrix}{f_{1} = {\frac{1}{N}{\sum\limits_{i = 1}^{N}\left( {d_{l} - D_{1}} \right)^{2}}}} & (4) \\{f_{1} = {\frac{1}{N}{\sum\limits_{i = 1}^{N}d_{l}}}} & (5) \\{f_{3} = {\frac{1}{N}{\sum\limits_{i = 1}^{N}d_{1}^{2}}}} & (6)\end{matrix}$

where D_(i) is the desired dose for the i^(th) sampling point and N thenumber of sampling points. Objectives of the form f₂ and f₃ are used forthe surrounding normal tissue, where the dose has to be minimized. f₂ isa measure of the integral dose, whereas with f₃ high dose values arepenalized stronger than with f₂.

For external beam radiotherapy and intensity modulated beam radiotherapy(IMRT) due to the sparse matrix nature of the kernel matrix we have abenefit only if the sampling point density is such that the averagenumber of sampling points per beamlet is larger than the number ofbeamlets.

Furthermore it is to mention that since the kernel of an inverseplanning system is based on the calculation of anatomy related objectivefunctions, it is obvious that the same benefit is expected when thisnovel method is implemented as described in the Equations 1 to 3 above.

In order to optimize for the specific organ's the objective function andthe derivatives are build up only by terms for which d_(i) ^(OAR)>D_(c)^(OAR)m_(S), which is expressed by Θ(d_(i) ^(OAR)−D_(c) ^(OAR)m_(s)). Alarge fraction of the sampling points has dose values that are smallerthan D_(c) ^(OAR)m_(S). It is possible to avoid the calculation of afraction of these dose values using the Cauchy-Schwarz-Bunjakovskiinequality:

${\sum\limits_{i = 1}^{N_{s}}{d_{i}b_{i}}} \leq {\sqrt{\sum\limits_{i = 1}^{N}a_{i\;}^{2}}\sqrt{\sum\limits_{i = 1}^{N}b_{i}^{2}}}$

For the dose d_(i) of the i^(th) sampling point we have

${d_{i}{\sum\limits_{i = 1}^{N_{d}}{x_{i}^{2}{\overset{\sim}{d}}_{il}}}} \leq {\sqrt{\sum\limits_{i = 1}^{N_{d}}x_{i}^{4}}\sqrt{\sum\limits_{i = 1}^{N_{d\;}}{\overset{\sim}{d}}_{il}^{2}}}$

The quantity

$r = \sqrt{\sum\limits_{i = 1}^{N_{d}}x_{l}^{4}}$

is calculated only once in each iteration while the N_(s) constants

$p_{i} = \sqrt{{\sum\limits_{i = 1}^{N_{d}}{\overset{\sim}{d}}_{il}^{2}},{i = 1},\ldots \mspace{14mu},N_{s}}$

can be calculated and stored in a pre-processing step. If r·p_(i)<D_(C)^(OAR)m_(s) it follows d_(i)<D_(C) ^(OAR)m_(s), otherwise it isnecessary to calculate d_(i). Even if the estimation of the dose by theinequality is not very good it is possible to avoid the calculation ofthe dose values of a large fraction of the sampling points in thespecific organs using only one multiplication per point.

A better estimate can be obtained using the relation:

${d_{i} = {{\sum\limits_{l = 1}^{Nd}{x_{l}^{2}{\overset{\sim}{d}}_{il}}} \leq {{\sqrt{\sum\limits_{l = 1}^{{Nd}\; 1}x_{l}^{4}}\sqrt{\sum\limits_{l = 1}^{{Nd}\; 1}{\overset{\sim}{d}}_{il}^{2}}} + {\sqrt{\sum\limits_{l = {{{Nd}\; 1} + 1}}^{Nd}x_{l}^{4}}\sqrt{\sum\limits_{l = {{{Nd}\; 1} + 1}}^{Nd}{\overset{\sim}{d}}_{il}^{2}}}} \leq {\sqrt{\sum\limits_{l = 1}^{Nd}x_{l}^{4}}\sqrt{\sum\limits_{l = 1}^{Nd}{\overset{\sim}{d}}_{il}^{2}}}}},{i = 1},\ldots \mspace{14mu},{Ns}$

where N_(d) is divided into 2 approximately equal terms N₁, N₂, i.e.N_(d)=N₁+N₂.

The two terms

$\sqrt{\sum\limits_{l = 1}^{{Nd}\; 1}x_{l}^{4}}\mspace{14mu} {and}\mspace{14mu} \sqrt{\overset{Nd}{\sum\limits_{l = {{{Nd}\; 1} + 1}}}x_{l}^{4}}$

can be calculated in each iteration once while the 2N_(s) terms

$\sqrt{\sum\limits_{l = 1}^{{Nd}\; 1}{\overset{\sim}{d}}_{{il}\;}^{2}}\mspace{14mu} {and}\mspace{14mu} \sqrt{\sum\limits_{l = {{{Nd}\; 1} + 1}}^{Nd}{\overset{\sim}{d}}_{il}^{2}}$

are calculate once and stored in a pre processing step before theoptimization.

This method can be extended to objective functions of the type

${f_{H}(x)} = {\frac{1}{N}{\sum\limits_{i = 1}^{N}{{\Theta \left( {{d_{i}(x)} - D_{H}} \right)}\left( {{d_{i}(x)} - D_{H}} \right)^{\alpha}}}}$

For a=2 we obtain the quadratic type of objectives, for a=1 theLessard-Pouliot objectives and for a=0 the DVH based objectives. Usingthe inequality we can avoid the calculation of the dose of a fraction ofthe sampling points if r·p_(i)<D_(H). For the rectum, bladder and normaltissue this approach avoids the calculation of a significant fraction ofthe dose values.

The main idea of the proposed dose optimization speed-up method is thatthe objective functions presented here and commonly used inbrachytherapy and their derivatives can be calculated without thecalculation of the individual dose values or a fraction of them. Forobjectives of the type given by equation (4)-(6) the number ofoperations for the calculation of the objective function values andtheir derivatives is independent of the number of sampling paints, if weignore the pre processing time to calculate and store {right arrow over(s)}^(S) ^(T) {tilde over (K)}^(S) and D^(S)={right arrow over (s)}^(S)^(T) {tilde over (K)}^(S).

This method allows us to increase the number of sampling points in theanatomical portion up to a few thousands, improving thus the accuracywithout any loss of optimization speed. In comparison to standard doseoptimization the new method is faster the more sampling points we have.The speed-up is also significant for implants with a small number ofsource dwell positions, where the term N_(s) ² is much smaller thanN_(s).N_(d) required previously for the calculation of the dose values.

A multi-objective optimization with 100 sources using up to 100solutions and up to 5000 sampling points in the anatomical portion witha 2 GHz PC is possible in less than 10 s. The storage for N_(s) samplingpoints and N_(d) sources can be reduced by a factor of approximately2N_(s)/N_(d).

A speed-up is expected not only for deterministic algorithms but alsofor stochastic algorithms such as genetic algorithms or simulatedannealing.

For high dose limits objectives an estimation of the dose value usingonly two multiplications and one addition per sampling point avoids thenecessity of the calculation of the dose value of a fraction of thesampling points.

FIG. 1 shows in very schematic form various elements of a known devicefor implanting an energy emitting source, e.g. radioactive seeds into aprostate gland. A patient 1 is shown lying in lithotomy position on atable 2. Fixedly connected to the table 2 is a housing 3. Housing 3comprises a drive means 4 to move rod 4 a stepwise. A template 5 isconnected or mounted to the table 2, which template is provided (notshown) with a plurality of guiding holes through which holes hollowneedles 9, 10 can be positioned relative to the patient. By means of aholder 6 a transrectal imaging probe 7 is fixedly connected to said rod4 a, which is moveable in a direction towards and from the patient bymeans of the drive means 4. The imaging probe 7 can be an ultrasoundprobe.

A needle 9 is used for fixing the prostate gland 11 in position relativeto the template 5. A number of needles 10 is fixed into position throughthe template 5 in the prostate gland 11. The template 5 determines therelative positions of the needles 10 in two dimensions. The needles 10are open at their distal ends and are sealed of by a plug ofbio-compatible, preferably bio-absorbable wax. In said housing 3 a seedloading unit 8 is present.

A well-known therapy planning module 12 a is provided for determiningthe number and relative positions of seeds in each needle forimplantation in the prostate gland 11. Such therapy planning module 12 ausually comprises a computer programmed with a therapy planning program.The therapy planning module 12 a is connected to the seed loading unit 8through a control device 12 for controlling the number of seeds for eachneedle. Control device 12 may be a separate device or may be anintegrated part either of the seed 1 oading unit 8 or of the therapyplanning module 12 a or may be embodied in the software of the therapyplanning module 12 a or of the seed loading unit 8.

The known device shown in FIG. 1 operates as follows. A patient 1 isunder spinal or general anesthesia and lies on the operating table 2 inlithotomy position. The (ultrasound) imaging probe 7 is introduced intothe rectum and the probe is connected via signal line 7 a with a wellknown image screen, where an image may be seen of the inside of thepatient in particular of the prostate gland 11 as seen from the point ofview of the imaging probe 7. The template 5 is attached to the drivemeans 4, thereby insuring the correlation of the ultrasound imagegeometry and the template 5. The prostate gland 11 is fixed relative tothe template 5 and the drive means 4 and the imaging probe 7 by means ofone or more needles 9, 10. Subsequently further needles 10 areintroduced in the body and the prostate gland under ultrasound guidanceone by one.

Moving the imaging probe with the drive means 4 longitudinally withinthe rectum controls the needle depths of each needle 10. After allneedles 10 have been placed, their positions relative to the prostategland 11 are determined in at least one of several known ways. In aknown way the therapy planning module 12 a determines how the needles 10are to be placed in the prostate and how many radioactive seeds are tobe placed in what order in each of the needles 10. The information aboutthe desired placement of the radioactive seeds in the needles 10 is usedto control the seed loading unit 8.

According to the invention said therapy treatment planning modulegenerates at least one treatment plan as it is provided with athree-dimensional imaging algorithm and a three-dimensional imagesegmentation algorithm for the specific organs within said anatomicalportion, the needles and the tubes for converting the image dataobtained with said imaging means into a three-dimensional image of theanatomical portion, whereby by using at least one single ormulti-objective anatomy based genetic optimization algorithm forpre-planning or virtual simulation purposes said processing means arearranged to determine in real time the optimal number and position of atleast one of said hollow needles, the position of said energy emittingsource within each hollow needle as well as the dwell times of saidenergy emitting source at each position; whereas for post planningpurposes said processing means are arranged to determine based onthree-dimensional image information in real time the real needlepositions and the dwell times of said energy emitting source for eachposition.

In FIG. 2 a template is disclosed for use in a real time radiationtreatment planning system according to the invention. Especially thetemplate 20 is detachable from a template frame 25, which frame isconnected with the stepper means for displacing the imaging means asdescribed in connection with FIG. 1.

According to the invention template 20 has a grid configuration withneedle holes 22 at an intermediate distance of 3.5 mm seen in diagonaldirection. In another embodiment the template has a grid configurationwith needle holes at an intermediate distance of 2.5 mm seen inorthogonal direction.

In a specific embodiment the template is a motorized template withoutholes and the needles are guided with a guiding tube, whereas theguiding tube can be positioned in each position of the virtual templategrid. As this embodiment does not use holes, the absence of a grid doesnot limit the positioning of the needles in relation to the template andthe anatomical portion to be treated. In fact with a template withoutholes the grid configuration is only limited to the diameter of theneedles used.

A more specific embodiment of the template is disclosed in FIG. 2, wheresaid template 20 is detachable from the frame 25. Frame 25 is connectedwith the stepper means as described above. For a good connection andorientation of the frame 25 and template 20 in relation to the device ofFIG. 1 the frame 25 is provided with alignment pins 27 which corporatewith corresponding openings (not shown) in the device of FIG. 1.

The template 20 has a saddle shaped body 20 a, which fits with the frame25 as shown in FIG. 2. For alignment purposes the template 20 isprovided with notches 11 a-11 b which cooperate with corresponding holes26 a-26 b present in the circumference of frame 25.

It is an another object of the invention to describe the catheters orneedles inserted in the body through which the HDR source is travellingwith their real geometrical dimensions. As a direct result of this it isa next object of the invention that sampling points for dose evaluation,which are lying inside the needles or catheters are excluded. This willcontribute to the reduction of the number of sampling points in theanatomical portion compared with other conventional methods and to theincrease of the speed.

As shown in FIG. 3 a catheter or hollow needle is defined by catheterdescribing points. These points are connected with cylinders and at eachcatheter describing lengths and diameters. The set of catheterscylinders and spheres are used to describe the geometry of a catheterthat may be either metallic linear or plastic and curved.

For generating each treatment plan the processing means of the radiationtreatment planning system according to the invention are arranged togenerate a set of multiple sampling points using said three-dimensionalimaging algorithm and said three-dimensional image segmentationalgorithm and to calculate the optimal radiation dose distribution foreach of said sampling points by using a gradient-based algorithm.

The quality of the results depends on the distribution of the samplingpoints as generated. According to the invention the dose distributioninside the anatomical portion (PTV), critical structures, such asspecific delicate organs (OAR) and the surface of the anatomical portionis estimated from the dose of a small number of points (samplingpoints).

As shown in FIGS. 4 and 5 the generated sampling points are distributedon the contours and on the triangulated surface of the anatomicalportion to be treated. For the contour based method no points are onboth ends of the anatomical portion. Therefore a large part of thesurface is undefined for the optimization algorithm and the resultingisodose is bounded only by the contours of the anatomical portion.

Sampling points in the volume are generated from low discrepancysequences or quasi-random distributed sampling points.

It is a another objective of the invention that in contrast topseudo-random distributed sampling points voids and clustering areavoided. Monte-Carlo generated quantities convergence much more rapidlythan a conventional pseudo-random sequence. Sampling points insidecatheters are excluded. This reduces the influence of very large dosevalues of sampling points that occasionally are produced very close tothe source dwell positions. Statistical values obtained from thesampling points are calculated therefore with a higher accuracy.

According to the invention two treatment planning steps are performed:Pre-planning or inverse planning: Given the geometry of the anatomicalportion (PTV) to be treated, the specific organs (OAR) near or withinsaid anatomical portion, a template and its position the optimal numberand position of needles, the dwell positions and the dwell times of theenergy emitting source are determined, so that the resulting dosedistribution satisfies various criteria such as coverage of theanatomical portion with the prescription dose, avoidance of dose valuesabove some critical values in the specific organs, etc.

Postplanning: Given the geometry of the anatomical portion (PTV) and thespecific organs (OAR) and a given number and position of needles and theposition of the energy emitting source in each needle the dwell times ofthe energy emitting source at each position are determined, so that theresulting dose distribution satisfies various criteria such as coverageof the anatomical portion with the prescription dose, avoidance of dosevalues above some critical values in the specific organs, etc.

Template based Inverse Planning

FIG. 6 defines the template and catheter characteristics. The planningsoftware is run on a personal computer or laptop computer and allows thesetting of certain objectives/boundaries/parameters prior to thegeneration of a treatment plan. The displacement step of the sourcewithin a needle is set at 5.0 mm in the afterloader parameters, since a2.5 mm value produces a large number of sources and the optimizationusing the algorithm according to the invention may take more time andsince then 512 MB RAM are recommended. By pressing the buttonAuto-activation the dialog of FIG. 7 appears, which may contain otherorgans (or VOIS i.e. Volumes Of Interest)

This dialog is used for the auto-activation algorithm. The anatomicalorgan (PTV) and the specific organs (OAR) to be protected against toomuch radiation exposure are listed as wels as the minimum distance ofthe source dwell positions from the corresponding VOI in mm. It is usedto select only source dwell positions that are at a distance to acorresponding VOI larger than a specified value.

VOIs for which the corresponding button is pressed only will beconsidered. In this example the rectum is ignored since it is outsidethe anatomical portion (PTV).

Now the program moves all catheters/needles, which inside the anatomicalportion taking the specific organs and the geometry of the anatomicalportion into account. The user has now to take only a subset of thesecatheters. In principle this will be done automatically by using anoptimization methods which are flexible and robust. This will be thetrue inverse planning.

In FIG. 8 the Source Parameter dialog for setting the prescription doseis disclosed. This dialog is used to define the source strength oractivity and the prescription dose. These parameters have to be suppliedfor the use of the optimization algorithms. The source is characterizedby its strength in units of U or as activity in units of GBq or Ci. Theprescription dose is specified in cGy.

Subsequently the Inverse planning dialog of FIG. 9 and FIG. 10 TemplateView and Loading appears. FIG. 10 shows the grid of the template, thecatheters and VOIs at various distances from the template. The selectedcatheters are shown in dark. The catheters which can be selected inlight gray. By moving the z-slide a plane parallel is moved along thenormal to the template and at some given distance from the templatedefined by distance z. In the Template View the intersection of the VOIswith that plane is shown in the anatomy window.

By selecting one of the buttons shown in FIG. 11 (and FIG. 10) thecatheter density can be selected.

By selecting with the mouse cursor one or more of the catheters orneedles (except those which are light gray) the selected catheter can beswitched on or off. So the user can select the catheters he wants to useduring treatment planning. For example a set of catheters on theperiphery and an additional set of catheters inside the anatomicalportion can be selected. It is preferred to limit the number of selectedcatheters or needles to 15-20 in order to limit the number ofcalculations to be performed.

Subsequently the dialog of FIG. 12 Geometry and sampling appears whichcontains information about the sampling points, the number of sourcedwell positions and the number of catheters. Subsequently a optimizationmodus has to be selected and the dialog of FIG. 13 appears. Here theoptimization method can be selected e.g. a deterministic optimizationmethod.

After selecting Set Optimization Options in FIG. 13 the dialog of FIG.14 appears. The specific organs (OAR=Organs At Risk) which are to beconsidered during the radiation treatment planning have to be selected.In this case the rectum is located outside the anatomical region to betreated and it can be ignored. However in this example prostate canceris to be treated and therefore the Urethra button is selected as theurethra is present inside the prostate. The critical dose value to whichthe specific organ may be exposed to as fraction of the prescriptiondose. In this case it is decided by the medical personnel that theurethra does not to receive more than 50% of the prescription dose.Therefore the factor 1.50 is entered.

With this dialog a single or multi-objective optimization can beperformed. The multi-objective optimization is selected. After theoptimization 20 solutions are presented to the user (medic personnel) ofthe treatment planning system. The treatment planning system Platodeveloped and commercialized by the present applicant Nuc1etron B.V. orother systems use a single set of importance factors which is notrecommended because there cannot be such a single set of importancefactors for all cases. There is not a single solution but in principleinfinite solutions.

The treatment p1anning system according to the invention tries toproduce a representative set of multiple treatment planning solutions.The deterministic method is the most simple approach. Recommended is ofcourse the evolutionary algorithm which is more flexible and producesmuch more solutions out of which the best for each case can be found. Itis not always possible to have similar results. Even if only prostatecases are considered.

After initialization the treatment planning system of the inventioncalculates the volumes of the anatomical portion and the specificorgans, it generates the sampling points and look-up tables are filled.In this case the optimization algorithm repeats 20 times with 20different sets of importance factors.

After the optimization step the Decision button has to be selected inorder to select˜solution and see the results. The dialog of FIG. 15appears. When the button Show results of all solutions is pressed thedialog of FIG. 16 appears. By moving the slider to the DVH values aremade visible, as these values are used in the decision making of thefinal treatment plan.

By selecting the column the DVH (1.500) urethra, the values in thatcolumn are sorted in descending order. See FIG. 18.

In FIG. 18 the best radiation dose coverage of the anatomical portion(PTV) is in this example 92.13%, while 10.785% of the urethra receives adose value above 1.5 times the prescription dose. If the medic personnelwant a dose exposure of the urethra below 1% (in FIG. 18 0.77%), thenthe best coverage for the anatomical portion (prostate) is 86.115%. Bypressing the Histogram button the distributions e.g. are displayed (FIG.19).

The deterministic algorithms use a mean on the dose normalization forthe surface of the anatomical portion and it is therefore not asflexible as the evolutionary algorithms. But the examples still show thedifferences between the treatment planning solutions obtained withdifferent importance factors/boundary conditions, which can be quitelarge. So one method would be to consider first the specific organs(OARs), then the dose coverage of the anatomical portion (PTV) andfinally the dose in the surrounding tissue. Whatever the preferences ofthe planner are the algorithm according to the invention generates allpossible solutions and the planner can select which treatment solutionis the best solution.

In the event that it is decided that 1% of urethra may to receive morethan the critical dose value, then the treatment solution no. 15 isselected in FIG. 20 (see also FIGS. 16 and 18). When the solution 15 inthe list is selected the Accept single solution button is to be pressedand for seeing the isodose distributions the Iso-Dose of SelectedSolution button is to be pressed in FIG. 20.

By selecting the 3D button in FIG. 21 the isodose values are marked,which are to be displayed (here the isodose for 1× the prescription and2× the prescription). Subsequently two 3D isodose distributions will bedisplayed.

Post Implant Optimization

Post Implant Optimization assumes that the source dwell positions aregiven. This is in principle what Nucletrons PLATO systems calls inverseplanning. After activating the Post Implant Optimization the treatmentplanning system loads the VOIS and catheters and the Autoactivationdialog of FIG. 22 appears.

After pressing the OK button the Source parameters dialog of FIG. 23 isdisplayed. After selecting the source parameters and pressing on OK thesystem directly continues with the optimization step of FIG. 13.Analogue to the pre-planning step the deterministic optimizationalgorithm can be selected. The steps of generating multiple treatmentsolutions are then the same as with the pre-planning step.

1-32. (canceled)
 33. A real-time radiation treatment planning system foruse in effecting radiation therapy of a pre-selected anatomical portionof an animal body, comprising: an imaging device configured to generateimage data corresponding to the anatomical portion; a radiation deliverydevice configured to deliver at least one energy emitting source throughat least one of a plurality of hollow needles inserted into theanatomical portion; a processor configured to generate, in real-time, aradiation treatment plan for effecting the radiation therapy, whereinthe processor is configured to: convert the image data into athree-dimensional image of the anatomical portion; designate a set ofmultiple sampling points in the three-dimensional image; calculateradiation dose distributions for the sampling points, based on at leastone of an actual number and position of at least one of the hollowneedles, an actual position of the at least one energy emitting sourcewithin the at least one of the hollow needles, or an actual dwell timeof the at least one energy emitting source at the position; calculate anobjective function indicative of an overall dose distribution in theanatomical portion, based on the radiation dose distributions of thesample points, wherein a number of operations for calculating theobjective function is independent of the number of sampling points; anddetermine the radiation treatment plan based on the objective function.34. The real-time radiation treatment planning system of claim 33,wherein the objective function is optimized using a gradient-basedmethod.
 35. The real-time radiation treatment planning system of claim33, wherein the processor is configured to determine, in real-time, theactual number and position of at least one of the hollow needles, theactual position of the at least one energy emitting source within the atleast one of the hollow needles, and the actual dwell time of the atleast one energy emitting source at the position, based on thethree-dimensional image and using an image segmentation algorithm. 36.The real-time radiation treatment planning system of claim 33, whereinthe processor configured to determine the radiation treatment plandetermines at least one of a virtual number, a virtual position, avirtual direction, or a virtual placement of the at least one of aplurality of hollow needles.
 37. The real-time radiation treatmentplanning system of claim 33, wherein the processor configured todetermine the radiation treatment plan determines a virtual amount ofradiation dose to be emitted.
 38. The real-time radiation treatmentplanning system of claim 33 further comprises a template configured toinsert at least one of the plurality of hollow needles inserted into theanatomical portion.
 39. The real-time radiation treatment planningsystem of claim 33, wherein the template includes a grid configurationhaving needle holes at a substantially constant intermediate distancealong a diagonal or an orthogonal direction.
 40. The real-time radiationtreatment planning system of claim 33, wherein the template includes aguiding tube for guiding the plurality of hollow needles, wherein theguiding tube is configured to be positioned as generated in thetreatment plan.
 41. The real-time radiation treatment planning system ofclaim 33, wherein the processor is further configured to calculate agradient function of the objective function, wherein a number ofoperations for calculating the gradient function is independent of thenumber of sampling points.
 42. The real-time radiation treatmentplanning system of claim 41, wherein the processor is configured toadjust the radiation treatment plan based on the gradient function ofthe objective function, such that the adjusted radiation treatment planyields an updated value of the objective function smaller than aprevious value.
 43. The real-time radiation treatment planning system ofclaim 33, wherein the sampling points are generated using a quasi-randomnumber generator.
 44. The real-time radiation treatment planning systemof claim 36, wherein the processor is further configured to: display thevirtual position of each of the plurality of hollow needles asdetermined in the treatment plan in the three-dimensional image of theanatomical portion on a screen; and guide the at least one of theplurality of hollow needles into the anatomical portion and to thevirtual position.
 45. A method for generating a radiation treatment planfor use in effecting radiation therapy of a selected anatomical portionof an animal body of a patient, comprising the steps of: generatingimage data corresponding to the anatomical portion; inserting at leastone of a plurality of hollow needles into the anatomical portion;converting the image data into a three-dimensional image of theanatomical portion; designating a set of multiple sampling points in thethree-dimensional image; calculating radiation dose distributions forthe sampling points, based on at least one of an actual number andposition of at least one of the hollow needles, an actual position ofthe at least one energy emitting source within the at least one of thehollow needles, or an actual dwell time of the at least one energyemitting source at the position; calculate an objective functionindicative of an overall dose distribution in the anatomical portion,based on the radiation dose distributions, wherein a number ofoperations for calculating the objective function is independent of thenumber of sampling points; and determining, in real-time, a radiationtreatment plan based on the objective function.
 46. The method of claim45, wherein the objective function is optimized using a gradient-basedmethod.
 47. The method of claim 45 further comprises determining, inreal-time, the actual number and position of at least one of the hollowneedles, the actual position of the at least one energy emitting sourcewithin the at least one of the hollow needles, and the actual dwell timeof the at least one energy emitting source at the position, based on thethree-dimensional image, or using an image segmentation algorithm. 48.The method of claim 45, wherein determining the radiation treatment plancomprises determining at least one of a virtual number, a virtualposition, a virtual direction, or a virtual placement of the at leastone of a plurality of hollow needles.
 49. The method of claim 45,wherein determining the radiation treatment plan comprises determining avirtual amount of radiation dose to be emitted.
 50. The method of claim45 further comprises calculating a gradient function of the objectivefunction, wherein a number of operations for calculating the gradientfunction is independent of the number of sampling points.
 51. The methodof claim 50 further comprises adjusting the radiation treatment planbased on the gradient function, such that the adjusted radiationtreatment plan yields an updated value of the objective function smallerthan a previous value.
 52. The method of claim 45 further comprises:displaying the virtual position of each of the plurality of hollowneedles as determined in the treatment plan in the three-dimensionalimage of the anatomical portion on a screen; and guiding the at leastone of the plurality of hollow needles into the anatomical portion andto the virtual position.